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Abstract 

Correlation effects of an electron gas in an external potential are derived using an 
Effective Action functional method. Corrections beyond the random phase approxima- 
tion (RPA) are naturally incorporated by this method. The Effective Action functional 
is made to depend explicitly on two-point correlation functions. The calculation is 
carried out at imaginary time. For a homogeneous electron gas, we calculate the effect 
of exchange on the ring diagrams at zero temperature and show how to include some 
of the ladder diagrams. Our results agree well with known numerical calculations. 
We conclude by showing that this method is in fact a variant of the time dependent 
density functional method and suggest that it is suitable to be applied to the study of 
correlation effects in the non-homogeneous case. 
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1 Introduction 



Correlation effects continue to be a popular subject due to the importance of these interac- 
tions in physical processes in general. Electron systems form an important subset of these 
systems. The nonlinearity that appears in this type of problem complicates the calculation 
of any physical properties of the system that strongly depend on the correlation between the 
individual particles. In this paper we are mainly concerned with a particular problem, that 
of calculating the ground state equilibrium energy of a large Fermi gas in an external poten- 
tial. At zero temperature the random phase approximation (RPA) is a well known method 
which addresses this problem. Here we present a way to improve on it, using an effective 
action method, which includes higher order interactions and which also applies at finite tem- 
perature and is applicable to non- homogeneous Fermi systems. Currently, such systems are 
frequently treated using the density functional method (DFT). PQ,[2] In DFT, the density 
of the system plays a central role and the energy is a functional of it. Both the method we 
present here and DFT are extensions of the Thomas-Fermi method. In DFT the density 
is found by self-consistently solving a one-particle Schrodinger equation. Correlation effects 
are taken into account by using the homogeneous result locally. The kinetic energy term 
is taken to be that of a non-interacting system having the same density as the real system 
under consideration. In spite of all these approximations, density functional methods have 
proved to give better values for binding energies than Hartree-Fock in atoms, molecules and 
solids in general. The method works best for nearly homogeneous systems; surface effects are 
handled poorly by this method. Attempts to include a gradient of the density in the energy 
functional gave mixed results. Similarly, systems with magnetic impurities do not give con- 
sistently good results when using DFT. The overall success of the method is therefore not 
well understood. Here we propose a functional approach in the spirit of DFT, that might 
offer advantages when treating inhomogeneity. The problem becomes harder to manage 
numerically in real cases. The method involves solving integro-differential equations which 
we solve in the homogeneous case. In the non-homogeneous case we suggest a possible 
strategy for a solution. 

Functional methods have given interesting results in many different areas of physics. 
0) |3] 0- |3 [l| They allow better treatment of non-perturbative effects and give a co- 
herent treatment of zero and non-zero temperature properties. A method such as the one 
presented here can be extended to study non-equilibrium properties of conducting ferromag- 
netic systems . jS] DFT is itself a functional method; Argaman and Makov gave a simple 
introduction to DFT based on this language. 0, (TO] Recently Kotani [H],[I2] proposed the 
use of correlation expressions based on the RPA approximation in ah initio calculations. His 
results were promising and consistent especially for magnetic elements. The method proved 
to be capable of giving better results than DFT in some cases. However, it was concluded 
that correlation expressions beyond the RPA are needed to give consistently better answers 
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than DFT. In this paper, we show how to get energy expressions which are more accurate 
than those obtained from RPA. Such expressions might be useful in numerical calculations 
such as the one carried out by Kotani. 11, We show below how to compute the free energy of 
the electron system. The system is assumed to be in an external static potential V(x). The 
result that we get is very general and could be applied to many different systems. It will be 
seen that the method provides a very compact expression for the energy that goes beyond 
the RPA method in a natural way. As an example, we apply our result for the energy to the 
homogeneous case in three dimensional space. We explicitly calculate the effect of taking 
exchange into account in the ring diagrams. Our results agree well with results obtained 
using quantum Monte Carlo methods. The treatment is done at finite temperature; JHj 
only at the end do we take the zero-temperature limit. Such results can be very useful for 
mesoscopic systems where correlation effects are important. 

The paper is laid out as follows: In section 2 we introduce the thermodynamic potential 
Q for a system subject to external sources. Q therefore will be a functional of these sources. 
The functional Q is the generator of connected Green's functions at finite temperature. 
Usually, one point source functions are used. Here, we will instead introduce two-point 
functions . The introduction of two-point external sources enables us to take into account 
the higher order corrections to the stationary phase approximation of the partition function 
in a simpler way. By including merely two diagrams (expressed in terms of the variables 
which are conjugate to the external sources) in the Effective Action, we are able to obtain 
in a very compact way the contribution to the energy of all the second order diagrams 
(and beyond) of the Coulomb interaction. This is the crucial advantage of this version of 
the Effective Action. The diagrams which are included here are two-particle irreducible, in 
contrast to the one-particle irreducible diagrams that we get when only one-point functions 
are used in the Effective Action. However, before bringing in the two-point correlation 
functions, we introduce a new field to replace the quartic Coulomb interaction by the well 
known Stratanovich-Hubbard transformation, so that all the interactions become local. This 
new field is simply the Hartree potential as in the RPA method. ^1] We explicitly introduce 
a term that describes correlation of the Hartree field at two locations. The sources can also 
be taken to be instantaneous so we can get a time-translation invariant solution. 

In section 3, we calculate the finite temperature Effective Action T. This action is 
obtained by a Legendre transformation from Q. Therefore, we introduce new variables 
conjugate to the sources. At the tree level, the T functional is simply the classical action 
of the system and is related to the free energy of the system. For homogeneous systems, 
T is known as the effective potential. In our case T will be a functional of three variables, 
the Hartree potential tp(r), the correlation function C(x,y) of the Hartree potential, and 
the Green's function p a /3(x,y) of the Fermi field. If the external sources are taken to be 
instantaneous, the points x = (ri,^) and y = (r 2 ,r^), correspond to the same temperature, 
T\ = r 2 , and the function p Q , j a(x,y) simply becomes the density matrix of the system, a 
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and (3 are spin indices. Even though it is possible to express the free energy in terms of the 
density matrix, it makes the manipulations more cumbersome. 

In section 4, we solve for <p(x) and C(x,y) in terms of p a p(x,y). We are able to do this 
because in reality p a p(x 7 y) is the only independent function of the problem. In this case we 
get an expression for the Effective Action V in terms of the density only. Therefore getting 
the expression for V in terms of p a p{x, y) is the essential result of this work from which other 
important calculations can be made. 

In section 5, we treat the zero temperature case. Because of the finite temperature 
approach, we have a set of diagrams, called 'anomalous', in the sense that they should 
not be present at zero temperature. The contribution of these diagrams vanishes at zero 
temperature. We go beyond RPA to include exchange effects on the ring diagrams. We also 
show that some of the ladder terms are included in our approximation. 

In section 6, we continue the treatment started in the previous section and calculate in 
detail the contribution of the ring diagrams if they include exchange. 

In section 7, we reexamine the non-homogeneous problem. We show how this method is 
related to the DFT formalism. In fact it is shown that a statement like the Hohenberg-Kohn 
theorem is trivially realized in this formalism. Similarly the question of v-representability 
can be given an affirmative answer within perturbation theory. 

Finally, the last section is the conclusion. 

2 The Thermodynamic Potential fl 

In this section, we obtain an expression for the thermodynamical functional Q to 'order' 
H 2 . [13] This functional is the logarithm of the partition functional Z of the system in the 
presence of external sources. The Z functional is written as an integral over all possible 
allowable paths of the different fields with a weighting factor that depends on the value 
of the action along the given path. The functional Q is simply the finite temperature 
generating functional of connected Green's functions of the fields involved. However because 
the Coulomb interaction, quartic in the Fermion field, is difficult to integrate, we introduce 
an auxiliary field by means of the well known Hubbard-Stratanovich transformation so that 
we can avoid the quartic interaction and its nonlocal behavior. In other words, the procedure 
consists of transforming the problem into one in which the electrons interact locally with 
a Hartree-type potential. In this treatment, exchange energy terms will show up directly 
because of the anti-commutation of the Fermi fields, which is built into the calculation as 
will be described below. The correlation terms are due to the quantum fluctuations around 
the Hartree potential. 

In the following we show how terms to order h 2 in Q can be obtained; we are including 
effects of second order in the full Hartree field of the theory. Therefore the equation of 
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motion satisfied by p(x, y) includes corrections beyond the Hartree-Fock approximation and 
involves second order exchange effects. 

A non-relativistic interacting electron gas in an external potential V(x) has the following 
Hamiltonian in the second quantized form: 



H n = / d 3 x 



-U>i (x)V 2 * a (x) + V(x)^ a (x)* a (x) 



/' 

+ l - J rf 3 ^^ft (l)$ t (2/) ^ (!/)$a(l) . (1) 

Here, we use units such that h = m = e = k B = l and j3 = ^. is a two-component 

temperature-dependent electron field Heisenberg operator, a and v are spin indices, i.e., 
a — 1 for spin up and a = —1 for spin down. Summation is implicit for repeated indices. 
In the following x will always mean a 3-D space vector. The system is constrained by the 
condition 

J d 3 x ¥ a (x)V a (x) = N (2) 

where N is the electron number operator which is constant. Therefore the density operator is 
simply ^^(x)^/ a (x) . Because of this constraint, we prefer to work instead with the following 
Hamiltonian 

H = H - fiN (3) 

where /i is a Lagrangian multiplier. The electron operator ^(x) can be expanded in terms 
of a complete orthonormal set of one-particle functions <fik( x ), so we write 

= ^2a k ip k (x), 

k 

k 

dk and a k are annihilation and creation operators. The subscript k includes momentum and 
spin. The ground state and the excited states can be represented as Slater determinants 
formed by the wave functions <Pk(x). In our case, it will be advantageous to take these 
functions to be self-consistent Hartree eigenfunctions. In the homogeneous case they reduce 
to plane waves. Since we are going to use a path integral formulation, we give the Lagrangian 
associated with the Hamiltonian H, 
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L 



d 3 x 



d 1 



if i 

— / d 3 xd 3 y— — -W(x)W(y)'$(y)'9(x). 
2 J \x-y\ 



(5) 



For a system at thermal equilibrium, we have to replace the time component t in the 
Hamiltonian by —it. The partition function Z is as usual defined as a sum over all possible 
states $. In the absence of external sources, we have 



Z[0] = £($|e-* H |$) 

= N J DVD&exp (-S[^,^ f ]) 



(6) 



where N is a normalization constant and the integration measure of the Fermi fields is only 
defined for fields that satisfy [To] : 



*(0,z) = -V(/3,x). 



(7) 



<S[\I / , vf"] i s the action of the system, 



drL 



dr / d 3 x < w(x 



\ fl 3 f l 

- / dr / d 3 xd 3 y — — -^ t (x)^ t (y)*(y)^(a;). 
2 Jo i 1^ 2/| 



(8) 



The functional integral for the partition function becomes 



Z[0] = N J DVD&expl- J drj 
+~ J dr dr y J d 3 y 



x 



^(z) 



F - 2/1 



<Kt-t„) 



(9) 



In writing the last term we made use of the fact that the Fermi fields satisfy a Grassmann 
algebra. An infinite self-energy term has been dropped from the above expression. Such a 
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term cancels at the end and has no effect on the evaluation of the eigenfunctions or related 
physical quantities. For notational convenience we set 



y) = (J^ - - 1* + s ( x - y) ( 10 ) 



and 



A{x, V )= S -^>. (11) 
\x-y\ 

Before proceeding further, we have to get rid of the quartic term in the Lagrangian. This 
is done by introducing an auxiliary boson field ip (x). We write 

Z[0] = N' J DVD&Dtp exp d A xd A y ^ a {x)G^{x, y)V a (y) 

~ J d A xd A y (f(x)A~ 1 (x,y)<f(y) + J d A x ^(x)^ (a;)* Q (a;)| . (12) 

The fourth component is the time component integrated over the proper range. It is easy to 
see that by using the following formula, 



/ 



dx e -^ xQx+bx = (det Q)"5 e ^ bQ ~ lb (13) 



and integrating over tp we get back the original expression for Z[0]. The prefactor N' is a 
new normalization constant. The operator A(x, y) is clearly invertible, 

A- 1 (x 1 y) = ~S(x-y)V 2 . (14) 

The non-local character of the Coulomb interaction has been removed by introducing the 
new bosonic field (p(x). It can be shown that <p(x) is the Hartree potential by using the new 
equations of motion derived from the new action of the problem. Now we couple the fields 
to local and non-local sources J(x),Q(x,y) and B(x,y). We first introduce the functional 
Q, a generator of connected Green's functions, through the normalized partition functional 
which is now a functional of the external sources. It is given by: 

Z[J,B,Q] = exp{-pSl[J,B,Q]} (15) 

such that 
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exp(-Q[J,B,Q]) J DVD&Dip exp (S[V, (p]) = 

J D*D&D(pexp^-S[V,V\<p] + J d 4 xd 4 y¥(x)Q(x,y)^(y) 

+ ^ d 4 xd 4 y<p(x)B(x,y)<p{y) + J d 4 x<p(x)J(x)} (16) 
where we have set (3 = 1. S[&, ip\ is simply the action of the transformed problem, 

S[*, <p] = J d 4 xd 4 y ^(x)G- 1 (x, y)tt(y) + l -<p{x)A-\x, y)<p{y) 

-J d 4 xip(x)^(x)^(x). (17) 



Now we define three new variables: 

6Sl[J,B,Q] . 
Jjjx) 



J=B=Q=0 = ((p(x)) = (p c (x) (18) 



\j=b=q=o = (*Ux)*p(v)) = Papix, y). (20) 



~^^yY ^ J=B=Q=0 = \ fa&My)) = \ \v<k x )<P<kv) + c ( x > y)} ( 19 ) 

6n[J,B,Q] 
5Q(x,y) 

(p c (x) is the expectation value of the field (p(x) in the ground state. C(x, y) is a correlation 
function of the field ip(x). p a p(x,y) is the Green function of the Fermi field. p a p{x,x) is 
therefore the density of the system. Here and below, the "time" ordering operator is not 
written explicitly. Therefore terms like p(x, y) are defined by setting t x — r y = + . Note that 
C(x, y) measures the departure from quasi-independence due to the correlation between the 
values of the potential at two different locations. The expectation value of the Fermi field 
is zero. 

In the rest of this section, we obtain an explicit expression for Q[J, B, Q] to 'order' h 2 . In 
the next section we solve for <p c (x), C(x, y) and p(x, y) in terms of J(x), B(x, y) and Q(x, y) 
to get an expression for the Effective Action T which is a functional of the new variables only. 
First we determine fi. We expand the exponent in the above integral around \I/ = ^ = 
and ip = <po where ipo is the configuration of ip that extremizes the action S. Therefore, we 
have 

(A' 1 + B)<p = -J. (21) 
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It is understood from the above that there is an integration over space and time on the 
L.H.S. of this expression. We choose from now on not to write integrals explicitly unless 
there might be some confusion. Now we expand around <po, so we write 



Viold) = V(new) + <A)- 



(22) 



Then assuming that the main contribution to the integral comes from the saddle point, we 
get the following expression for the partition functional: 



exp { - W = exp(- fl 5W) det|5 " + Q|det|A " + B|4 E- 



x exp 

where we have defined S and E to be 
S = exp 



detfG- 1 ] dctlA- 1 ]-^ 
J D^D^D^ E j^e 2 (<^ f *) 2 +...| 



-tf^G" 1 + Q)V - -ifiA- 1 + B)tp 



1 



-^(G- 1 + Q)tt - -ifiA' 1 + B)<p 



det M 



S = y DVD&Dip exp 
In the above, we have used the fact that 

da cfot e 

for Grassmann numbers a and at, and 

J dx e~^ xMx = (detM)- 1 * 
for c-numbers x. The argument of the first term on the right is 



and 



Now using the fact that 



1 1 



G~ = G~ - v^o- 



(23) 

(24) 
(25) 



(26) 



(27) 



(28) 



(29) 
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det^ = e TrlnA , (30) 
we can get an explicit expression for fl[J, B, Q\. 

Q[J, B, Q] = S(<fi ) - Tr ln[l + G + GQ] + \ Tr ln[l + A B] 

— \DmD^D^(^m) 2 - + ... (31) 
2 J S 

After finding Q, we now calculate the Effective Action T. We have to solve for the external 
sources in terms of the physical variables p, C and (p c . This we do in the next section where 
we find the Effective Action at finite temperature. However, this treatment applies to zero 
temperatures too. In both cases we get an expansion in h . However h is only used in 
bookkeeping the diagrams included in the expansion of Q. 

3 The Effective Action T 

The Effective Action is a functional of <p c (x), C(x,y), and p(x,y), which is obtained by a 
triple Legendre transformation from Q[J, B, Q] 

T[(p c ,C,p] = Q[J, B,Q] - J d 4 x ip c (x)J(x) - ^ J d[xy] ^ c (x)B(x,y)tp c (y) 

~\ J d \ x y\ C(x,y)B(x,y) - J d[xy]p(x,y)Q(x,y) , (32) 

where we have written d [xyz...] for d A xd 4 yd A z... to simplify the notation. It is understood 
that the integration over time is carried out for r in the range [0, /3] for the nonzero temper- 
ature case and over all time for the zero temperature case. 
It is easily verified that : 



5T[^C,p] 
8ip c (x) 



-J(x)- J d 4 y B(x,y)<p c (y) (33) 



5C(x,y) 2 

and 



5r[<f c ,C,p\ 1 

■^B{x,y), (34) 



5p(x,y) 
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When we turn off the external sources, the above equations give the values of f c ,C and p 
that minimize the Effective Action. (p c (x) and C(x,y) are really dependent variables since 
they depend on the auxiliary field <p(x), so we should be able in principle to express them 
in terms of p(x, y) which gives us the density. 

Now we have to express J, B and Q in terms of ip c , C and p. First we note that when 
h = 0, we have 

sn . . 

Jj = (36) 

Hence, we can write 

fo = ¥>c + <p (37) 
where tp is of order H. Similarly, we write 

S(<p )=S(<p c )+hS 1 . (38) 
Now we seek an expression for F in the form 

r = T + hT 1 + h 2 T 2 . (39) 

It then follows that 

r = -vv4~V- (40) 

Now, we find approximate expressions for B and Q in terms of (p c , C and p. Using Eq. (j31j) 
and Eq.()35p. we get the following relation: 

p- 1 G= 1 + GQ -eGp c + 0{h). (41) 
The expression for B is more involved. Again using Eq. (J35|) and treating (p as a functional 
of J and B, we get: 

BC = (l-2e P (p e ) (i- 2 ^) -A- X C. (42) 



Inserting back all these expressions in T, we get for V 



i 



T 1 = -Tr In (p^G) + ^Tr ln(C~M) + ^Tr (A^C) - Tr (pG^ 1 ). (43) 
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The terms of order H 2 are also straightforward but more involved. The steps are similar 
to those in the case of the Effective Action with one-point sources only. ^7] Actually, it 
can be shown that the next terms in T 2 are the sum of two-particle irreducible diagrams of 
the theory. [3] , [S] Here, we sum only the first two diagrams in this series expansion. Fig. ^ 
shows the first four diagrams in this expansion. The first two diagrams enable us to include 
first and second order exchange effects. Diagram e is not part of T 2 since it is two-particle 
reducible. However it is one-particle irreducible and it is part of the usual Effective Action, 



Figure 1: Graphs a, b, c and d are part of the expansion terms in T. Graph e is however 
reducible and does not show up in the expansion. The solid line represents the propagator 
p(x,y), the dashed line represents the propagator C(x,y). 

Hence the reducible graphs do not appear in the expansion. It is the term E that appeared 
in Eq. (|23|) that is responsible for excluding such a graph from the expansion. T 2 and higher 
order terms represent vacuum graphs where the full propagators are p(x,y) and C(x,y) of 
the Fermi field and the boson field, respectively. Therefore to order h 2 , the Effective Action 
has the following expression 



[14j 




(44) 



S{<p c ) + -Tr In AC' 1 + A~ X C 
Tr \np- 1 G + G- 1 ^ c )p + T 2 [C,p} 



(45) 
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where 



d_ 



-V 2 - n + V(x) - ip c {x) 



(46) 



The term jj, J d 3 x p(x, x) cancels the term /iN in the Hamiltonian H, and therefore it can be 
handled without difficulty. However, we have to deal with the \x that appears in the term 
Tr In G. This can be removed by a convenient choice of the path of integration in the 
complex u;-plane. The case \i = occurs when there are no charges present. We are dealing 
with negatively charged electrons bound by a static potential V(x). Therefore bound states 
appear, with negative energies bounded from below for any physically sensible V(x). In 
Fig. El we show how by going from the C\ path, which corresponds to nonzero /jl, to the Co 
path we pick up contributions from the bound states between the two paths, since energies 
of bound states are poles of the propagator G^ 1 in the complex cu plane. Note that the 
eigenvalues Uk of these states do not take into account correlations due the Coulomb field. 
Because of the logarithmic operator, there is a cut along the real positive axis. By going 
back to real time and then Fourier transforming the time dependence, we get 



where tr applies only to spatial variables. The uVs are the eigenvalues of the equation 



The single particle functions ipk(x) of the potential V(x) are not supposed to be used alone 
as the starting point for any numerical calculations since they do not take account of the 
repulsion between electrons. However, it is clearly advantageous to separate the effect of the 
external potential V(x). This separation also appears in other treatments such as DFT. 

The ground state energy E g is given in terms of the Effective Action per unit time. Using 
the above results, and the fact that the time independent values of p(x, y), C(x, y) and <p c (x) 
are defined by setting t x = r y = 0, the expression for E g is 




(47) 




(48) 
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Figure 2: Path of integration used to obtain Eq. (J47j) . 



E g J dt = \J d l x v] tPcWA-^yfyciy) + l -Tr lnC A' 1 
-- / d[xy] A~\x,y)C(y,x) - u k I dt 

-Tr \npG~ 1 + J d A x (J- - l -W 2 + V(x) - e ip c {x)\ p{x,y)\*=l Ty 

+^e 4 J d[xyuv] p(x, y)p{y, u)p(u, v)p(v, x)C(x, u)C(y, v) 

+7>e 2 J d[xy] C(x,y)p(x,y)p(y,x) (49) 

where we have dropped terms of order higher than e 4 in T 2 for simplicity. This is the full 
expression for the ground state energy of the system. Here G _1 has p set to zero. To be 
able to use this expression for E g , more approximations are required. Since f c is a dependent 
field, in principle we should be able to solve for ip c and C in terms of p only. This we do 
next, but only approximately to keep the calculations manageable. 
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4 The Effective Action as a Functional of p(x, y) 



In this section, we find an expression for T solely in terms of the Green's function p(x,y). 
This is done by finding the expressions for (p c (x) and C(x,y), that minimize T, in terms of 
the density p(x,y). For (p c (x) we get the following expression: 



73„. p{y,y) 



ip c (x) = -e / d 6 yf^r. (50) 

f — y\ 



Similarly, minimizing T with respect to p(x,y), we get 
d 1 



z) 



__ _ 1 V 2 + V(x) - e tp c (x 



p(x,z)+e 2 j d 4 y p(x,y)C(x,y)p(y, z) 
e 4 / d[yuv]p(v,u)p(u,x)p(y, v)p(z,y)C(u,y)C(v,x) . (51) 



Finally, minimizing T with respect to C, we get 



5{x,z) = j d A y[A 1 {x,y) - e 2 p{x,y)p{y,x)] C{y, z) 

-e 4 / d[yvu}C(u,v)C(y,z)p(x,u)p(u,y)p(y,v)p(v,x) . (52) 



To be able to solve for C(x,y) we have to linearize the theory. In the following, we keep 
terms only to 'order' e 4 . We find that the correlation C(x, y) is given by 



C(x,z) = A(x,z) + e 2 j d[xix 2 } A(x,Xi)D(xi,x 2 )A(x2, z) 

+e 4 / d[xix 2 x' l x' 2 } A(x 1 x' l )D(x' l ,x' 2 )A(x' 2 ,xi)D(xi,X2)A(x2,z) 

where we have set 

D(x l ,x 2 ) = p(xi,x 2 )p(x 2 ,xi) . (54) 

The above equation for C(x,z) can be represented diagrammatically as shown in Fig. El 
The first term is the bare Coulomb potential. The second and the third are the direct and 
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Figure 3: The expansion of the propagator C(x, z) to order e 4 



exchange term with one bubble. So far the Fermion propagators are the true propagators in 
this expansion. In the following we will ignore the corrections with two bubbles and higher. 

Clearly C(x,y) is the full screened Coulomb potential if all orders of e 2 are included. 
Similarly, to order e 4 , the equation satisfied by p(x, z) is 

JL _ 1 v 2 + V{x) + e 2 / dy^l] p(x, z) = 6(x, z) 
_ ut x z j | x y | _ 

-e 2 / dyA(x,y)p(x,y)p(y,z)-e 4 / c%aax 2 ] p(x,y)p(y, z)A(x,x 1 )D(x 1 ,x 2 )A(x 2 ,y) 



+e 4 y d[yxix 2 ] p(x 2 ,a;i)p(xi ! ?/) / o(x,X2)p(y,2;)yl(xi,x)A(a; 2 ,y) . (55) 

This is the equation satisfied by the Fermi Green's functions .^1] The term on the R.H.S. 
of order e 2 is an exchange term or Fock term due to the statistics of the electrons. The 
last two terms of order e 4 take into account collisions and are equivalent to the usual Born 
approximation for two-particle Green's functions in scattering theory, Fig. 0] 

The most likely way to solve these integro-differential equations is by iteration. Using 
the above equations, we can write an explicit expression for the energy in terms of the full 
propagator p(x, z) of the theory. However, such an expression suffers from the problem of 
over-counting some of the states of the system. This is mainly due to the fact that e 2 is not 
a good expansion parameter. Hence physical arguments are used to discard some terms at 
this level of approximation and instead include them at higher orders. 

An expression for the electron propagator in terms of the free propagator is given in 
Fig. El • The corresponding analytical expression can be easily written following usual rules 
jTJJ. Here we keep only the first three terms. The second term vanishes in the homogeneous 
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Figure 4: Born approximation for two-particle Green's function. 



Figure 5: Approximate solution to the one-particle Green's function. The propagators on 
the right correspond to the free theory. 

case, the case we are mainly interested in this paper. We expect the expression for the 
energy, which has not been linearized, to be a good one if we believe that a stationary phase 
approximation is viable. Without including T2, the expression obtained for the energy is 
correct with p(x, y) the Hartree propagator and ip c {x) the Hartree potential. This indicates 
that a stationary phase approximation is possible and that the higher order correction T2 
should be a small perturbation to the Hartree solution. The linearization of the problem 
must take account of this. Hence, we can justify the validity of the expansion of Y in h, 

r = r + hTi + h 2 v 2 + ... (56) 
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by the result we get in the end. To second order in h, Eq.(|53|) becomes 
dyA-\x,y) C(y,z) = 5(x,z) + he 2 J dyD(x,y)A(y, z) 

+he 4 J d[yuv] p(x, u)p(u, y)p{y, v)p(v, x)A(u, v)A(y, z) 
-2ft 2 e 4 / d[yuv]p(x,u)p(u,v)p(v,y)p(y,x)A(u,v)A(y,zl57) 



where now p(x, y) is the non- interacting electron propagator. We have used only a first 
order approximation to the true propagator to find this equation for C(x,y). Using this 
equation, we get an expression for the energy to order h 2 with proper symmetry factors for 
the diagrams involved in the expansion. Setting h = 1, we obtain for the zero temperature 
limit, 
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TE g [p] 



dx dy 



p H (x,x)p H (y,y) 



\x-y\ 

+ J dx ^-^V 2 p H (x,y)^j \ x=y 
1 



+ e / dx V(x)ph(x, x) 



UJ k >fl 



- e 4 / dxidx 2 dxdy p H (xi,x)p H (x,x 2 )pH(x2,y)pH(y,x 1 )A(x,y)A(x2,x 1 ) 



+ - Tr In 

2 



S(x,z)+e 2 J dyA(z,y)D H (x,y) + 



dydx x dx 2 p H {xi, x)p H (x, x 2 )p H (x 2 , y)pii(v, x 1 )A(x 2 , Xt)A(y, z) - 



2e 4 J dydxidx 2 p H {x, x 1 )p H (x 1 , x 2 )p H {x 2 , y)pn(y, x)A(x 2 , x{)A{y, z) 



—Tr In 



8(x, z) -e 2 dy p H (y, y)p H (x, z)A(x, y) - e 2 dy A(x, y)p H {x, y)pii(y, * 



-e 4 J dydx x dx 2 p H (x,y)p H (y,z)p H (x 1 ,x 2 )p H (x 2 ,x 1 )A(x 1 ,x)A(x 2 ,y) 
+e 4 / dydx x dx 2 p H (x 2 ,x 1 )p H (x 1 ,y)p H (x,x 2 )p H (y, z)A(x 1 ,x)A(x 2 ,y) 



e 2 dy p H (y,y)p H (x,z)A(y,x)\ x=z + e 2 / dxdy A(x,y)p H (x,y)p H (y, z) 



— e 4 y dxdydx\dx 2 b{x — z) A(x, xi)A(y, x 2 )pn(y, z)ph(x 2 , xi)pn{xi, y)pH(x, x 2 ) 

+e 4 J dxdydxidx 2 5(x - z)p H (x, xi)p H (x u x 2 )p H (x 2 , y)p B {y, x)A(x 2 , xt)A(y, z)(5 

T is an interval of time. The trace, Tr, acts on x and z. So now we have obtained an 
expression for the energy in terms of the Green's functions p(x, z) only. This is a major 
result of the work. As we mentioned above, we solve for p(x, y) to order e 4 by iteration of 
Eqs. (|51|52j) . The expression for the energy can be solely written in terms of p(x, x). However 
as we stated earlier this is not advantageous and it is better to keep on working in terms of 
p(x,y). 

The classical Coulomb term and the interaction with the external potential appear nat- 
urally in this approximation. They can be easily separated from the full expression for the 
energy. Another important point is that a gradient of the density also appears naturally 
within the above expressions. By expanding the first logarithm, we immediately obtain the 



18 



Hartree-Fock exchange term, i.e., 



-e 2 / dxdy 



p H (x,y)p H (y,x) 



S{r x -T y ). 



(59) 



'exch 



x-y 



The usual expression for the exchange energy follows trivially from Eq. fl58|) upon using the 
commutation relation of the Fermi field. The infinity arising from this cancels the self- 
interaction term that was dropped in Eq.®. There is a similar exchange term that comes 
from the second logarithmic term . This term is canceled by another similar term outside 
the logarithm. The usual RPA term is contained in the first logarithm. The new, extra 
term in the same logarithm provides among other things exchange corrections to the ring 
diagrams. The higher order exchange diagrams are also included in this approximation. We 
will say more on this when we treat the homogeneous case in the following section. 



5 The Homogeneous Electron Gas at Zero Tempera- 



In this section we apply the main result of section 4, the expression for the energy Eq. ()58)l . to 
the homogeneous case at zero temperature. The literature on this problem of the calculation 
of the correlation energy is huge. The most complete treatment, so far as we are aware, 
was given by Bishop and Luhrmann who used a linked cluster type of expansion which 
was restricted to zero temperature effects. ^H] We assume that there is a background of 
positive charge with density equal to the average density of the electron gas, i.e., the system 
is neutral. Since the system is homogeneous, the final expression for the energy will be given 
in momentum space. The Green's function that will be used in the following is the solution 
to the first iteration of the nonlinear equation satisfied by p(x, y). Since there is no external 
potential, the input Green's function is that of a free electron gas. The energy expression is 
given in imaginary time, hence we use the following expression for the free Green's function 



In this section and the next, we use this definition of the electron Green's function which 
differs by a factor of i from the one used previously. In terms of po(x, y) the ground state 
energy is given by the following (omitting the subscript for simplicity) 



ture 




(60) 
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T E gs 



where 



x=y 



J d 4 x (^iV 2 p(x,y)*j 

e 4 f e 2 f 

+— / dx dz 5{x — z) Ai(x, z) + — / dx dy p(x,y)p(y,x)A(x,y) 

~Tr M In [8(x,z) - e 2 H(x, z) + e 4 M(x,z) - 2e 4 IT(x, z)] 



+7V(:r,z) In [^(a;, -2) _ e 2 IC(x, z) — e 4 S(x, 
+ y dxdz5(x — z) <y— TC(x, z) — — .A4(x, z) + e 4 E(x, 2) + e 4 n(x, 2;) j- (61) 



,M(x, z) = J dydxidx2p(y,X2)p(x2,x)p(x,xi)p(xi,y)A(xx,X2)A(y,z), (62) 

ft(x,z) = J dyA{z,y)D(x,y), (63) 

n(x, *) = / dy d Xl dx 2 p(y, x)p(x, x l)p ^ x 2 )p(x 2 , y)A^ x 2 )A(y, z) , (64) 

tC(x,z) = J dyA(x,y)p(x,y)p{y,z), (65) 

S(x, 2) = Si (a;, 2) - E 2 (x, 2), (66) 

£1(2,2:) = J dydx 1 dx 2 p( y x,y)p(y,z)p( y x 1 ,X2)p(x2,x 1 )A(x,x 1 )A(x 2 ,y), (67) 

E 2 (x,z)= / dydxidx2p(x,x 2 )p(x2,xi)p(xi,y)p(y,z)A(xi,x)A(x2,y). (68) 



The above equation, Eq. (j61|) . is simply Eq. (j58|) taking into account of the fact that p(x,x) 
gets canceled by an equal and opposite charge. Because of the translational invariance of 
the system, we can write a simple expression for the correlation energy in the momentum 
representation, 
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p d p ' dco ' 

E c = |eY J Yl ( 27r )3 2tt A(kPi ~ P ^ A< ^ 2 ~ Vz)p{V\)p{V2)p{Vz)p{V\ ~ P2 + Pa) 

+\v J {In [1 + e 2 H{p) + e 4 M(p) - 2e 4 U(p)} - e 2 H{p) - e*M(p) + 2e 4 n(j9)} 

-V j { In [1 - e 2 IC(p) - e 4 S(p)] + e 2 IC(p) + e 4 £(p) } . (69) 



Since the ring diagrams are part of E c , we can test Eq. fl69|) to see if our expansion gives the 
correct leading result. First we show how this equation follows from Eq.()61|). We start by 
letting 

H{x, z) = - J dy D{x, y)A{y, z) (70) 
or, since the system is homogeneous, we can write instead 

H(x-z) = J dyD(x-y)A(y-z). (71) 

Next, we rewrite D(x — y) and A(x — z) in terms of their corresponding Fourier transforms, 
i.e., 

/dk 
—-- exp[ik ■ (x - y))D(k) (72) 



and 



My - z) = j j^y 4 exp[zg • (y - z)}A(q). (73) 



Hence, by convolution we have, 

Hip) = (A*D) (p), (74) 

however since 

d(p) = j j0r 4 p(p + q)p(q), (75) 
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it follows then that, 



Hip) = A(p) J j^y 4 p( P + q)p(q). (76) 



After using the expressions for the propagators, we end up with the following expression for 

nip) 

with q — (u>,q). Integrating first over u) in the complex plane, we find that 

H(v) = -2e 2 ^ I {W+ qr ~ jgj ^ - ~ " lg± (78) 

where we added a factor of 2 to account for spin. This term H(p) will prove to be all that 
is needed to reproduce the known RPA result. All other terms that appear in the energy 
expression are new additions to the correlation. Before showing this explicitly , we give the 
expressions for the remaining terms JC(p), Ai(p), n(p), £i(p), £ 2 (p) and 

JC(P)=P(P) J j^pip + q)A(q) (79) 
M(p) = A{p) j ^-^ p (k 1 )p(k 1 +p)p(k 2 )p(k2-p)A(k 1 - k 2+P ), (80) 

n(p) = A(p) | ^ J ^ P (k 1 )p(k 1 -p) P (k 1 -p- q )A( q ), (si) 

^P^^p(Pi)p(P2)p(P2+Pi -p)-4(p 2 -p)^(p-pi), (82) 
S 2 (P) = p(p) / 7^7^^(Pi)p(p-Pi)p(Pi+P2)p(p 2 ) (83) 



and 



S(p)=p(p) | ( ^ ) > (fei)p(fe 2 )p(fc 2 + fc 1 -p) ^(fcx-p^-ACp-feMCp-fcx)]. (84) 
The Fourier transform A{p) is given as usual by 
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Figure 6: Second order exchange diagram in Eq. EH1 



A(p) 



Att 



(85) 



The terms in Eq. flfi9j) have the following meaning after we expand the ln-terms. The first 
term represents the second order exchange term, Fig. |3 The term 7i(p) is responsible for 
generating the ring diagrams, Fig. [7| The term Ai (p) is responsible for generating some of 
the ring diagrams with exchange. The cross terms give the remaining ring diagrams with 
exchange, Fig.lHl The H(p) term generates ring diagrams with a self-energy insertion in each 
ring, fig. E3 Hence the terms in the first In should contribute the most. Some of the terms 
that appear upon expansion of the second In term are shown in Fig. E3 The term K,(p) 
generates the so called 'anomalous' diagrams, Fig. ^2[13j The term £i(p) generates ladder 
terms like those in Figs. ^] and ED ^{p) generates terms like those in Fig. El Cross terms 
of the last two terms are shown in Fig. ED Diagrams that appear in Fig. EDare important 
in a nonhomogeneous medium and must be accounted for since they no longer vanish. They 
arise whenever p(x, x) is not equal to the ion density in the energy expression, Eq.(jSBJ)- 

Normalizing the momentum with respect to kp, the Fermi momentum, the explicit ex- 
pression for the correlation energy per particle in Rydbergs (Ryd.) has the following form, 



N 




n(x,y)-M{x,y)+2U(x,y)} 




dx / dy y 2 {In [1 + K.(x, y) + E(x, y)} - K.(x, y) - E(x, y)} , (86) 
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Figure 7: Ring diagrams generated by TC(p) in Eq. (|69|). 



Figure 8: Diagrams that appear in the first In term in Eq.( 169)1 . 



Figure 9: Ring diagrams generated by II (p) in Eq. 
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Figure 10: Diagrams that appear in the second In term in Eq.( IfiHj) . 



Figure 11: A new infinite set of diagrams at nonzero temperature. 
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Figure 12: A term of order e 4 that is due to Ei. Both representations are equivalent. Initially 
two particles interact with each other with one of them going back to its initial state while 
the other one interacts with a third particle before both returning to their corresponding 
initial states. 



Figure 13: Some of the ladder diagrams due to E^ 
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Figure 14: Several terms that appear due to S 2 . 



Figure 15: A term that is both due to Si and E2. Here three particles are interacting 
pairwise with two of the three particles exchanging states. 



Figure 16: Coulomb interactions that must be taken into account in inhomogeneous media. 
These diagrams vanish in the homogeneous case. 
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where we have set E qi = \q + Pi\ and E q2 = \q — p 2 \. 
For convenience, we set 

_£ = I, + j 2 + / 3 



(87) 



where Ji, J 2 , and J 3 are the first, the second and the third integrals, respectively. The 
functions K,(x,y), H(x,y) and Yli(x,y) are given by, 



lC(x,y) 



4ar s g(y) 



n(—ix + y 2 — 1) 



= -e 2 A(x,y)D(x,y) 



or more explicitly, we have 
2ar 



«<*,«0 - ^h^(^ + V(i-|)(i + |))in(| 



|)2 + y 2 (1 + p 



"2y V 



)2 + y2(l-K)5 



— ( arctan f ^ 2 J~ — — arctan ^ - 



and 



Si(a;,y) = a / d 3 k 1 d 3 k 2 g(x,y,k 1 ,k 2 ) 



where 
g(x,y,k u k 2 



(89) 



(90) 



(91) 



i i x (eg - h) - eg - k l2 )) (ejh - g) - g(i - k 2 )) 



y- k^ \y- k 2 \ 2 (~ix + y 2 - l)(ix - k\ - k\ + E> - 1) 

9(x) is the step function. E', g(y) and a are given by 

E'= \h + k 2 -y\\ 



(93) 



f 4 



and 



(94) 



g(y) 



2y 



1 + 2/ 



(95) 



28 



The Gell-Mann-Brueckner term can be isolated from the full expression for the correlation, 
Eg. ( 186(1. by rewriting the second integral, I2, in the following way: 



where 



h 



ring . jexch , jseRing 



jring 



rexch 



3 1 

4"7T r 2 a 2 
3 1 
47r r 2 a 2 



■00 
00 



dx \ dy y 2 {In [1 + H(x, y)] - H(x, y)} 



dx I dy y < In 



1 + 



M(x,y) 
1 + H(x,y) 



-M(x,y) 



(96) 

(97) 
(98) 



rseRing 



is what is left of To, — /™ n9 — If*"*. The term J" n9 is indeed the full expression for 
the RPA term j^E^ PA . By taking the limit r s — > and y < 1, it reduces to the Gell-Mann- 
Brueckner result 



E? PA 
N 



0.0622 In r s -0. 142 + ... 



as we show next. Going back to the expression for 7i(p), Eg.( I7H|) . we have 



(99) 



jring 



with H(k,uj) defined to be 
<u(i \ o 2 47r f d3 P 

Now we set 



Z7T 



(2tt) 3 



Kp + ^) 2 - \f 



k 2 J (27T) 3 (j(^ + £)2_ l£2 



2e 2 4tt 



(2tt) 3 F 



(H 1 (h,u) + H 2 (k,u)) 



(100) 



(101) 



(102) 



Terms other than J"" 9 in Eg.f 186(1 provide corrections to the RPA-term. The first term in 
Eg.(IEEl), as we mentioned above, is the second order exchange term,[Tj 



E 2 r exch = -In 2 ((3) 

3 2vr 2SV ' 



0.04836ityd. (103) 



29 



To compare our expansion to others we calculate the term that is equivalent to the RPA 

calculation. We start by evaluating the integral Hi. After integrating the angular variables, 
we have 



Hi(u>, k) — 7r / — — In 



u 2 + (pk + \k 2 ) 2 
u 2 + {-pk + \k 2 ) 2 



(104) 



This integral can be easily performed over p. We get after that 



Hi(co,k) 



7T r i 



A; 3 [ 2 

-im 

2 



-fc 2 W 



In 



^2 + ^2(^ + ^)2 



u; 2 + A; 2 (A; F --)(A; F + -) 



+ k 2 (-k F + ^) 2 

, fc(A* + §)\ ( H-kp + l 

arctan — — arctan 



+ k F k 3 } . (105) 



The second integral H 2 in H(u>, k) is seen to be simply obtained from Hi(cu, k) by replacing 

— * — * 

k by —k. Therefore we have 



H(u,k) = 2Hi(u, k) 



(106) 



If we set y = and x = -, we get 



2ar s 

7T?/ 2 



^ + y(i-V )ln 

(2y) 3 



x 2 + V(l + |) 2 
x 2 + V(l-|) 2 



arctau I — 2 ^ + arctan ^ 



x 



x 



(107) 



where r s = Hence the ring diagrams' contribution to the energy is given by 

3N 



j^RPA _ 



Ana 2 i 



y 2 dy dx[\n[l + H(x,y)]-H(x,y)\, 



(108) 



where N is the total number of electrons. Now, we notice that if we make the following 
substitutions: 
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Figure 17: Ring diagrams with exchange. 



y = y' 

x = 2yx' (109) 

we recover exactly the same expression for the ring diagrams as obtained by Bishop and 

Luhrmann (BL) using a totally different expansion. ^H] This helps to validate our original 
expansion in h and later in iterating the equations of motion in terms of e 2 . We stress that 
this expansion is valid for both small and large r s . 

From the above analysis, we see explicitly that the method of Effective Action amounts 
to including another infinite set of diagrams besides the usual ring diagrams. One subset of 
the diagrams added is the ring diagrams that allow exchange in them, Fig. |H| On physical 
grounds these are expected to be the next important ones that must be summed up. It is 
also easily seen from the above that all second order diagrams are included in the expansion 
with the right symmetry factors. Once more, this shows that our original expansion in h 
is indeed meaningful. One last thing to note about the diagrams in Fig. El is that they 
include the "anomalous" diagrams which appear due to the finite temperature method. [21] 
From a nonzero temperature calculation, we were able to show that these diagrams give 
a zero contribution at zero temperature simply because of the Fermi statistics. From 
Eq. fjHUj) the contribution of these diagrams involves integrating In [1 + K(x, y)\ —)C(x, y) over 
all x in the complex plane. Since JC(x,y) has a simple pole, Eq.([88|). the above integrand 
ends up having poles of order two and higher, resulting in their zero contribution at zero 
temperature. However these diagrams become essential at non-zero temperature where the 
above constraint of Fermi statistics is no longer an issue. In the following section we calculate 
the contribution of the function Ai(x, y) to the correlation energy. 
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Figure 18: Comparison of the Hubbard approximation and the BL approximation to the 
Coulomb exchange interaction : (-= — i — — T ) vs. y . 

to x \k 1 -k 2 +y\ 2/ y 



Figure 19: Comparison of final results for the correlation energy (in Ryd.) of a homogeneous 
Fermi gas: 1. Ref . |23| . 2. Ref. |18j . 3. This paper 



Figure 20: Diagrams that appear when we include terms of order H in T. 
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6 The Inclusion of Second Order Exchange Effects 
in Ring Diagrams 

In this section, we give the contribution of diagrams like those shown in Fig. El to the corre- 
lation energy at zero temperature. Hence we go beyond RPA in this case. We will show that 
our method provides excellent agreement with fully numerical calculations. We also compare 
our results to those found by BL.[TH] where a coupled cluster formalism has been used to get 
the correlation energy. To get to these final results we had to make approximations along 
the way. We use two different approximations: the Hubbard approximation [22] and the 
Bishop-Luhrmann approximation. JB] We have found that the former applies well to high 
values of r s while the latter applies well to low values of r s . From Eq.(|98|). this amounts 
to finding I£ xch and the second order exchange term that has already been calculated by 
Onsager et al.^H] Now we show how to calculate Jf xc/ \ The calculation is straightforward 
but special care must be exercised when it comes to numerically evaluating the final result. 
First we let 

/oo /*oo 
dx / y 2 dy In 
■oo Jo 

where 



1 + 



M(x,y) 
1 + H(x,y) 



(110) 



A<(*,y) = ^l/^£*» 1 7 ^ (111) 

1 y> (2tt) 4 y* J r \ h -k 2 + yr (_ ix + M+|!z^ U + Mz|!zMl 

and where the region of integration T is given by: 

r = k 2 )h <l,k 2 < 1, \y+ h\ > 1, \y- k 2 \ > l} . (112) 

It is obvious from the above that we are faced with a daunting task of having to deal with 
a 9-D integral inside a logarithmic function which itself has to be integrated over normalized 
energy x and normalized momentum y. Most of the complications are related to the angle 
integrals which can not be separated. To be able to make some progress we have to make an 
approximation, i.e., we either assume that \ki + y\ ~ 1 and average over the angle between 
them, i.e., the Hubbard approximation or we use a more sophisticated approximation like 
the one proposed by Bishop and Luhrmann (BL). Hubbard's approximation amounts to the 
following: 

— -. 113 

\h-k 2 + y\ 2 V 2 + l 
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On the other hand, the BL-approximation is more complicated and the reader is referred 
to their appendix for a discussion of their approximation. Both approximations are 
almost equivalent for large values of momentum y, i.e. y > 2 in units of Fermi momentum. 
For values of y less than 2, both approximations are quite different ( see Fig. El )• The 
Hubbard approximation seems to apply well for moderately large values of r s while the BL- 
approximation applies well for low values of r s . In fact the BL-approximation incorporates 
the Fermi statistics of the interacting particles where it matters most, i.e. in high density 
situations, and this is one reason why we found it to be a good approximation to the exchange 
Coulomb interaction for low r s . The Hubbard approximation should be expected to be a 
good approximation for particles near the Fermi surface and where statistics are not an issue. 
This case should apply well to low density. Given either approximation we end up with an 
expression for Ai(x,p) of the form: 

lQa 2 r 2 s 1 ^ f x 2 + \{y 2 + 2k\ ■ y)(y 2 - 2k 2 ■ y) 



M(x,y) = -—f-f(y) / d*hd%- C v (114) 

( 27r ) V J r [x 2 + \(y 2 + 2k 1 -y) 2 j (^x 2 + \(y 2 -2k 2 -y) 2 j 

where the function f(y) is our approximation to the exchange term in Ad(x, y) and is assumed 
known. The integrations over k\ and k 2 are easily performed if we choose to write both 
vectors in cylindrical coordinates with y along the z-axis, i.e., 

h = (pi,Zi,6i), i = 1,2 . (115) 
Hence the integration over k\ becomes: 

/■27T p 

d?ki = I d6i I pi dpidzi, (H6) 



o Jr 



where 

r' = [y<2, -l <Zl <l-y, (l-( Zl + y) 2 )^ < Pl < (l-^ 2 )'} 
\j{y<2, l-y< Zl <l, 0<p 1 <(l-z 1 2 )"} 

\j{y>2 1 -1< Z1 <1, 0< Pl <(l- Zl 2 )^, (117) 

with an equivalent expression for the region of integration over k 2 . The integrations are now 
easily carried out. The expression that we get for Ai(x, y) is naturally expressed in terms of 
Y instead of x, so in the following x refers to We only quote the final expression here: 

cy 2 r 2 1 

M(x, y) = -16^- f(y)£( X , y), (118) 
ir z y^ 
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C(x,y) 




H(x,y) (H(x,y)-2) 




arctan 



2-V 



x 
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( X 2 + ( 2 + p f) 
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Given this expression for .M(x, y), we can easily write the full expression for the correlation 
due to exchange effects in ring diagrams. It is important to note that we allow for more 
than one exchange to take place in each ring diagram. Hence our original integral I, can be 
now evaluated numerically. So within the above approximation, we were able to reduce our 
calculations to a two-dimensional integral: 



jexch _^ 

2 47r (ar,) 2 



dx 



y 



In 



1 + 



M(x,y) 
1 + H(x,y) 



M(x,y) . 



(120) 



The expression for the function f(p) depends on which approximation we choose to use. 
In this formalism it is hard to decide on which one since both give very close answers to the 
second order exchange diagram. We have found that the Hubbard approximation in this 
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instance gives 0.041 Ryd. while the BL approximation gives 0.049 Ryd. Clearly the latter 
is closer to the true value of 0.0484 Ryd., but we believe that this is not enough to decide 
against using the Hubbard approximation. In fact at low r s we expect that the high energy 
part of the exchange energy between particles to dominate and this explains why we found 
the BL-approximations to apply well in this range. At high r s this is not true anymore and 
in this instance we expect the BL- approximation to overestimate this exchange of energy 
between the particles. This is in fact what we have found. This creates doubt that the 
BL-approximation is always better than the Hubbard approximation. As we stated above, 
the second order exchange term is not a good test since it is density independent. This 
explains why the Hubbard approximation gives better results for higher values of r s . Below 
we give some numerical answers to compare with those obtained using quantum Monte Carlo 
methods (QMC). [231 Our results are obtained using both approximations. For the second 
order exchange diagram we have used the true value of 0.0484 Ryd. 

Clearly our results, FigUHl agree well with previous results. For large values of r s , i.e., r s 
> 10, our results do not compare well to the QMC results. However this is expected since we 
need to include corrections due to the ladder diagrams which are believed to be important 
at low densities. This concludes our calculations. We have shown the effectiveness of an 
Effective Action approach to electronic systems. In the next section we shall examine the 
viability of using this formalism in the nonhomogeneous case and show how it is related to 
density functional theory. 

7 The Nonhomogeneous Electron System 

In this last section, we go back to the inhomogeneous case and show how our method gives 
some of the results known in DFT. £[] We start by briefly reviewing the basic ideas 
behind DFT. The Hamiltonian H is written as the sum of three terms: T, V and U. T is the 
kinetic energy term, V is the external potential term and U is the Coulomb energy term, so 
we have 



H 



T + V + U 



(121) 




(122) 




(123) 




(124) 
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If \l/ is the ground state, then the density n(r) is given by 

n(r) = {9f\ip\r)ip{r)\^f). (125) 

The first important fact to note is that V(r) is a unique functional of n(r), assuming that 
there is a unique ground state for the system. Now let the functionals F[n{r)} and Ev[n(r)} 
be defined as follows: 

F[n(r)] = {*\T + U\V), (126) 

E v [n(r)] = J V(r)n(r)d 3 r + F[n(r)}. (127) 

Then, if we assume that there is a one-to-one correspondence between n(r) and V(r), it can 
be shown that 

T^/U-O. (128) 

with 

J n(r) d 3 r = N . 

This result is called the Hohenberg-Kohn theorem. Next, the Coulomb energy is isolated 
from the functional F[n] by introducing a new functional G [n\. 

F[n\ = - [ drdr^^P-^p- + G[n] , (129) 
2 J \r — r I 

G[n] = \J d ? V r V^ ni(r,r')| r =r' + ^ ^ y~7[ dr rfr '' ( 130 ) 

C 2 (r, r') = n 2 (r, r'\ r, r') — n(r)n(r'). (131) 



and 



Here ni(r, r') is the one-particle density matrix and C^r, r') is a correlation function defined 
in terms of the one- and two-particle density matrices. DFT calculations are essentially 
centered around finding good approximations to the functional G[n\. This is usually done 
by postulating that there is a virtual system of free electrons in an external potential with 
exactly the same density as the interacting system. The energy is found by finding the 
eigenfunctions that correspond to this external potential self consistently. 
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The method that was presented makes use of concepts which are similar in many ways 
to the ideas expressed in DFT. However there are major differences. Our method clearly 
incorporates the Hohenberg-Kohn theorem. From Eq. (j35|) . we have the following 



Sp(x, x) 



-Q{x,x). (132) 



From Eq.(|49J). the term 



/ 



d 3 x V(x)p(x, x) , (133) 



can be separated from the energy E gs . Hence if we set Q(x, x) = , we immediately get 
the Hohenberg-Kohn result. In fact we now have an explicit expression for -EVb 7 -] within 
perturbation theory. This result also shows the one-to-one correspondence between density 
and external potential as long as there is a unique solution around (J,B,Q) = (0,0,0) for 
the defining equations, Eqs. (J33H33J) . 

We have solved these equations only approximately. We have made two important 
approximations. The first corresponds to the number of diagrams included in T from the 
outset. However from the calculations of the homogeneous case, we see that it will probably 
be the case that including only two diagrams in the inhomogeneous case will give good 
results. Including higher order corrections by taking into account diagrams like the ones 
in Fig. 1201 is nevertheless straightforward if more accurate results are desired. The second 
approximation we have made was in solving for p(x,y). We had to linearize the equations 
of motions to be able to solve for p(x, y) iteratively. Also from Eq. (joT?j) . we get the following 
equivalent result: 

VVc(r) = 4vm(r). (134) 



Clearly the method presented here is a generalization of the Thomas-Fermi method where 
n(r) in the last equation is the true density of the system. A very important difference from 
DFT is that within this method we have a systematic scheme for calculating the functional 
G[n], Eq. (j49|) . From the above analysis in the homogeneous case, it is obvious that the 
most important contribution to correlation at zero temperature comes from calculating the 
following term: 



A' 1 + e 2 ppX + e A XXppApp 



Tr In X = In det -r « ~r* ™ HH"HH {m) 



where 

X = A- X C . (136) 
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To calculate the determinant of the above operators, we need a basis of wavefunctions. In 
practical computations, it is more advantageous to use a finite basis set. If we choose as our 
basis set the wavefunctions {(j)i}f =1 with eigenvalues such that 

±\7 2 + V(x^J <f> i {x) = e i <l> i {x), (137) 

then the wavefunctions, if)j, of the full interacting system of N electrons can be written as a 
linear combination of the above wavefunctions: 

N 

ipj(x,t) = Y, a jM^)e~ ezt ■ (138) 

i=l 

The coefficients dji are found self consistently by solving Eq. (JH5J) . To solve for X(x,y), 
we solve Eq. (j52|) self-consistently. Clearly this method is computationally more intensive 
than DFT. However approximations like the static approximation used by Kotani JT] can 
simplify the calculations a lot. 

Finally, we would like to point out that this method applies equally well to excited states 
other than the ground state. The Hartree field of the excited state is found again by using 
Eq.tjHHl). If ipo(x) is the ground state solution then ipi(x) = <p c (x) + A(x) is the excited state 
solution iff 



8(p(x 
This implies that 



5T 

<p=<pi 



0. (139) 



A f , J .,A(y) = 0- (140) 



8<p(x)6<p(y) 



In the above, A (a;) is assumed small compared to tpo(x). This latter equation enables us to 
solve for A(x) and hence the corresponding p(x, y). We expect this method to apply equally 
well to atomic and molecular systems as we have shown to be the case for homogeneous 
systems. 



8 Conclusion 

We have used a functional method, the Effective Action method, to calculate correlation 
effects in electronic systems. A general expression for the energy was derived which includes 
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corrections beyond the RPA approximation. We have also shown how to apply this result to 
a homogeneous electron gas. Our results agree well with quantum Monte Carlo calculations. 
The inclusion of higher order corrections showed that the energy expression derived is not 
necessarily applicable only to high density cases. This method has many similarities with 
the principles on which Density Functional Theory is founded. The main result is that 
a systematic expression for the correlation energy can be written down. The use of such 
expressions in realistic calculations, with approximations like those used by Kotani 
should result in our being able to simplify the calculations, not only in the homogeneous 
case but also for the nonhomogeneous case. 
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r s QMC-1 



0.1 
1.0 
2.0 
5.0 
10.0 



0.121 
0.0902 

0.0563 
0.0373 



BL-2 



GRPA-3 



0.245 0.243 

0.123 0.121 

0.0930 0.0949 

0.0612 0.568 

0.0429 0.0419 



